/*------------------------------------------------------------------------
 *  Copyright (C) 2011  Luis M. de la Cruz Salas
 *
 *  This file is part of TUNA
 *
 *  TUNA is free software: you can redistribute it and/or modify
 *  it under the terms of the GNU General Public License as published by
 *  the Free Software Foundation, either version 3 of the License, or
 *  (at your option) any later version.
 *
 *  TUNA is distributed in the hope that it will be useful,
 *  but WITHOUT ANY WARRANTY; without even the implied warranty of
 *  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
 *  GNU General Public License for more details.
 *
 *  You should have received a copy of the GNU General Public License
 *  along with this program.  If not, see <http://www.gnu.org/licenses/>.
 ------------------------------------------------------------------------*/

#ifndef _GENERALMESH_H_
#define _GENERALMESH_H_

#include <vector>
#include "Tuna.hpp"

namespace Tuna {

/*! 
 ***************************************************************************
 * General class for meshes. 
 * This class is the root of the hierarchy for mesh classes. This should be
 * an "abstract" class from the point of view of the design. The functions 
 * \c calcPoints() and \c calcElements() should be implemented in the 
 * derived classes.
 ***************************************************************************
 *  \author  Luis M. de la Cruz [ Sat Jan 19 18:16:37 GMT 2010 ]
 *************************************************************************** 
 */
template<typename Tmesh>
class GeneralMesh {
    
public:
/*! Precision used for the numerical operations. */
  typedef typename Typeinfo<Tmesh>::prec_t prec_t;     
/*! Arrays for numbers*/
  typedef typename TunaArray<prec_t,1>::huge ScalarField1D;

  GeneralMesh() { 
    Dim = Typeinfo<Tmesh>::Dim; 
    xyz[0].resize(1); 
    xyz[1].resize(1); 
    xyz[2].resize(1); 
  }
  ~GeneralMesh() { }

/*!    
 * This funtion gets an returns a reference to the derived class,
 * using the Curiously Recursive Template Pattern (CRTP). 
 */
  Tmesh& asDerived() { return static_cast<Tmesh&>(*this); }    

/*! 
 *  Delegate responsabilities to derived classes. 
 */
  bool calcPoints()   { return asDerived().calcPoints(); }  

/*! 
 * Delegate responsabilities to derived classes. 
 */    
  bool calcElements() { return asDerived().calcElements(); }  
/*!
 *  Returns a std::vector<prec_t> with the coordinates distribution
 *  in the axis defined by the parameter \c a.
 *  \param a can be \c X, \c Y or \c Z  
 */
  ScalarField1D getCoord(Axis_t a) 
  {
    ScalarField1D aux( xyz[a].size() );
    int i = 0;
    for(viter = xyz[a].begin(); viter != xyz[a].end(); ++viter, ++i)
      aux(i) = *viter;
    return aux;
  }

  inline vector<double> getX() { return xyz[0]; }
  inline vector<double> getY() { return xyz[1]; }
  inline vector<double> getZ() { return xyz[2]; }
  
/*!
 *  Write the coordinates of the mesh to a file. The file is in ASCII
 *  and 1, 2, or 3 columns are written depending on the dimension
 *  1D, 2D or 3D respectively.
 *  \param filename is a string for the name of the output file.
 */
  void writeToFile(string filename) {
    std::ofstream file (filename.c_str());
    if (!file) 
      std::cout << "\n +-----+"
		<< "\n | TUNA::GeneralMesh : Can't open file \" " << filename << "\""
		<< "\n +-----+";   
    if( Dim == 1 ) 
      for(unsigned int i = 0; i < xyz[0].size(); ++i)
      	file << xyz[0][i] << "\t" << 0 << "\n";
    if ( Dim == 2 ) 
      for(unsigned int j = 0; j < xyz[1].size(); ++j) {
      	for(unsigned int i = 0; i < xyz[0].size(); ++i)
      	  file << xyz[0][i] << "\t" << xyz[1][j] << "\t" << 0 << "\n";
      	file << "\n";
      }

    if (Dim == 3 ) 
      for(unsigned int k = 0; k < xyz[2].size(); ++k) {
      	for(unsigned int j = 0; j < xyz[1].size(); ++j) {
      	  for(unsigned int i = 0; i < xyz[0].size(); ++i)
      	    file << xyz[0][i] << "\t" << xyz[1][j] << "\t" << xyz[2][k] << "\t" << 0 << "\n";
      	  file << "\n";
      	}
      	file << "\n";
      }
  }

protected:
  int Dim;               ///< Dimension of the mesh.
  vector<prec_t> xyz[3]; ///< Coordinate arrays.
  //  ScalarField1D xyz_coord[3];
  typename vector<prec_t>::const_iterator viter;

private:

};

} // Tuna namespace 


#endif //_GENERALMESH_H_



